関数定義ほか

単日データを処理するための関数を定義しておきます。

関数名 処理概要 必要なパッケージ
lagdiff(n) 単日データの前日差を求める dplyr
ma7(n) 単日データの7日間(1週間)移動平均を求める zoo
ma28(n) 単日データの28日間(4週間)移動平均を求める zoo
ma(n, k) 単日データのk日間移動平均を求める zoo
lagdiff <- function(n) {        # 前日差を求める関数
  n - dplyr::lag(n, default = 0L)
}

ma7 <- function(n) {            # 移動平均(7日)を求める関数
  zoo::rollmeanr(n, k = 7L, na.pad = TRUE)
}

ma28 <- function(n) {           # 移動平均(28日)を求める関数
  zoo::rollmeanr(n, k = 28L, na.pad = TRUE)
}

ma <- function(n, k = 7L) {     # 移動平均(k日)を求める関数
  zoo::rollmeanr(n, k = k, na.pad = TRUE)
}

 

日次集計のための関数を定義しておきます。この関数は第一引数にデータフレームを取り、データフレームが以下のような構成を取る個票データであることを前提としています。

df <- data.frame(
  date,   # Date型の日付データ
  key,    # グルーピングのためのデータ(例:都道府県、地方区分など)
  ...     # 順不同
)

ここで定義してる集計関数は任意の変数の水準ごとに日次集計を行います。なお、エラー処理の記述については省略しています。

daily_aggregate <- function(df, date, key) {
  date <- dplyr::enquo(date)   # NSE(Non-standard evaluation)処理
  key <- dplyr::enquo(key)     # NSE(Non-standard evaluation)処理

  df %>% 
    dplyr::group_by(!!date, !!key) %>%    # 任意の水準ごとの日次集計を行う
    dplyr::summarise(n = dplyr::n()) %>%  # n()は個数をカウントする関数
    dplyr::ungroup() %>%                  # 最後にungroupするのがポイント
    tidyr::complete(                      # 暗黙の欠損を補完する
      date = seq.Date(from = min(!!date), to = max(!!date), by = "day"),
      !!key, fill = list(n = 0L)          # 個票がない=陽性者ゼロ
    ) %>% 
    dplyr::group_by(!!key) %>%            # 上記の日次集計を元に各種計算を
    tidyr::nest() %>%                     # 行うには nest & unnest 処理を
    dplyr::mutate(                        # 利用する
      diff = purrr::map(data, ~ lagdiff(.$n)),   # 前日差
      cum = purrr::map(data, ~ cumsum(.$n)),     # 累計
      ma7 = purrr::map(data, ~ ma7(.$n)),        # 移動平均(7日)
      ma28 = purrr::map(data, ~ ma28(.$n))       # 移動平均(28日)
    ) %>% 
    tidyr::unnest(cols = c(data, diff, cum, ma7, ma28)) %>% 
    return()
}

この関数は tidyr::nest, tidyr::unnest を用いてグループごとにまとめた日次集計データを purrr パッケージを用いて処理してる点がポイントです。後述の全国集計と何が異なるのかは実際にコードを動かして確認してください。

その他、グラフで繰り返し使う文字列を変数として定義しておきます。

subtitle <- paste0("Generated @", lubridate::now())
caption <- "Data Source: covid19japan.com"

 

Import & Tidy

個票データの集計に限らず、データをインポートした際には目的に応じて各変量(変数)の変数型を変更します。特に水準ごとに層別処理を行いたい場合には因子型に変換しておくと便利です。また、結合したいデータと名称や体系を合わせておくこともポイントです。

 

都道府県データ

prefs <- "https://gist.githubusercontent.com/k-metrics/9f3fc18e042850ff24ad9676ac34764b/raw/f4ea87f429e1ca28627feff94b67c8b2432aee59/pref_utf8.csv" %>% 
  readr::read_csv() %>% 
  dplyr::mutate(
    # Googleの予測データと結合するためにコード体系を合わせる
    japan_prefecture_code = paste0("JP-", `コード`)
  ) %>% 
  dplyr::select(
    # Googleの予測データと結合するために名称を変更する
    japan_prefecture_code, prefecture_name = pref,
    # 日本語の変数名は扱いにくいので英語名に変更する
    pref = `都道府県`, region = `八地方区分`, pops = `推計人口`
  ) %>% 
  dplyr::mutate(
    # 水準ごとに表示させるために因子化する(あらかじめデータをコード順に
    # 並べておくことが因子化の際のポイントのひとつ)
    japan_prefecture_code = forcats::fct_inorder(japan_prefecture_code),
    pref = forcats::fct_inorder(pref),
    region = forcats::fct_inorder(region),
    pops = as.integer(pops)
  )

prefs

 

Covid19Japan個票データ

df <- "https://raw.githubusercontent.com/reustle/covid19japan-data/master/docs/patient_data/latest.json" %>% 
  jsonlite::fromJSON() %>% 
  dplyr::select(patientId, date = dateAnnounced, gender,
                pref = detectedPrefecture, patientStatus, knownCluster,
                confirmedPatient, ageBracket,
                deceasedDate, deceasedReportedDate) %>% 
  dplyr::filter(confirmedPatient == TRUE) %>% 
  dplyr::mutate(
    date = lubridate::as_date(date),
    gender = forcats::as_factor(gender),
    pref = stringr::str_to_lower(pref),
    patientStatus = forcats::as_factor(patientStatus),
    cluster = dplyr::if_else(!is.na(knownCluster), TRUE, FALSE),
    ageBracket = forcats::as_factor(ageBracket),
    deceasedDate = lubridate::as_date(deceasedDate),
    deceasedReportedDate = lubridate::as_date(deceasedReportedDate)
  ) %>% 
  # 都道府県データと結合
  dplyr::left_join(prefs, by = c("pref" = "prefecture_name")) %>% 
  dplyr::select(-pref) %>% 
  dplyr::rename(pref = pref.y) %>% 
  # 因子型の欠損値を水準化する
  dplyr::mutate(
    japan_prefecture_code = forcats::fct_explicit_na(japan_prefecture_code,
                                                     na_level = "JP-48"),
    pref = forcats::fct_explicit_na(pref, na_level = "空港検疫"),
    region = forcats::fct_explicit_na(region, na_level = "空港検疫"),
    gender = forcats::fct_explicit_na(gender, na_level = "非公表"),
    ageBracket = forcats::fct_explicit_na(ageBracket, na_level = "非公表"),
    patientStatus = forcats::fct_explicit_na(patientStatus,
                                             na_level = "Unknown")
  )

df

 
作成したデータフレームを確認します。意外と欠損値(n_missing項参照)が多いデータであることが分かります。

df %>% 
  skimr::skim()
Data summary
Name Piped data
Number of rows 394269
Number of columns 14
_______________________
Column type frequency:
character 2
Date 3
factor 6
logical 2
numeric 1
________________________
Group variables None

Variable type: character

skim_variable n_missing complete_rate min max empty n_unique whitespace
patientId 0 1.00 1 16 0 394269 0
knownCluster 391764 0.01 3 88 0 233 0

Variable type: Date

skim_variable n_missing complete_rate min max median n_unique
date 0 1 2020-01-15 2021-02-02 2020-12-20 371
deceasedDate 393889 0 2020-02-13 2020-11-19 2020-05-08 151
deceasedReportedDate 393939 0 2020-02-13 2020-10-17 2020-05-16 131

Variable type: factor

skim_variable n_missing complete_rate ordered n_unique top_counts
gender 0 1 FALSE 3 非公表: 284597, M: 61370, F: 48302
patientStatus 0 1 FALSE 9 Unk: 391735, Hos: 1261, Dec: 372, Hom: 315
ageBracket 0 1 FALSE 13 非公表: 284694, 20: 29433, 30: 19042, 40: 16089
japan_prefecture_code 0 1 FALSE 48 JP-: 100858, JP-: 44112, JP-: 41222, JP-: 25631
pref 0 1 FALSE 48 東京都: 100858, 大阪府: 44112, 神奈川: 41222, 埼玉県: 25631
region 0 1 FALSE 9 関東地: 203051, 近畿地: 77883, 中部地: 40024, 九州地: 34494

Variable type: logical

skim_variable n_missing complete_rate mean count
confirmedPatient 0 1 1.00 TRU: 394269
cluster 0 1 0.01 FAL: 391764, TRU: 2505

Variable type: numeric

skim_variable n_missing complete_rate mean sd p0 p25 p50 p75 p100 hist
pops 2174 0.99 7915.48 4232.27 560 5107 7537 13822 13822 ▆▅▆▇▇

 

Data Wrangling

 

陽性者判定者数

 

全国集計

日次集計を行います。集計は dplyr::group_bydplyr::summrise を用います。日付データがない日は集計として現れませんので、tidyr::complete でデータを補完しておきます。
その後、集計結果を元に前日差、累計、移動平均を計算します。

japan_daily <- df %>% 
  dplyr::group_by(date) %>% 
  dplyr::summarise(n = dplyr::n()) %>% 
  dplyr::ungroup() %>% 
  tidyr::complete(
    date = seq.Date(from = min(date), to = max(date), by = "day"),
    fill = list(n = 0L)
  ) %>% 
  dplyr::mutate(
    diff = lagdiff(n),   # 前日差
    cum = cumsum(n),     # 累計
    ma7 = ma7(n),        # 移動平均(7日)
    ma28 = ma28(n)       # 移動平均(28日)
  )

japan_daily
変数 変数の意味 変数型 備考
date 陽性判定者の報告日 Date 陽性判定日とは異なる
n 陽性判定者数 Integer
diff 同前日差 Integer 初日は陽性判定者数に同じ
cum 同累計 Integer
ma7 同移動平均(7日) Double
ma28 同移動平均(28日) Double 4週間(1ヶ月)

 

地方区分別集計

ここでは定義した daily_aggregate 関数を用いて集計を行っています。処理としては全国集計の場合と同じ考え方ですが、前述のように前日差や累計などを計算する際のデータの処理方法がことなっている点に注意してください。

region_daily <- df %>% 
  daily_aggregate(date, region)   # 上記の集計処理を関数化したもの

region_daily

 

都道府県別集計

地方区分集計と同様に集計します。

pref_daily <- df %>% 
  daily_aggregate(date, pref)

pref_daily

 

年代別集計

 

全国集計

ageBracket_daily <- df %>% 
  daily_aggregate(date, ageBracket)

ageBracket_daily

 
年代別の地方区分別集計ならびに都道府県別集計は割愛します。理由はVisualize章のグラフを見て頂ければ分かると思います。  

クラスター別集計

 

全国集計

cluster_daily <- df %>% 
  daily_aggregate(date, cluster)

cluster_daily

 
クラスター別の地方区分別集計ならびに都道府県別集計は割愛します。理由はVisualize章のグラフを見て頂ければ分かると思います。  

Visualize

作成した集計データを可視化します。

陽性判定者数

全国集計

japan_daily %>% 
  dplyr::filter(date == max(date))
subset <- japan_daily
title <- "【全国】陽性者数(単日)"
xlab <- ""
ylab <- ""
sec_scale <- 50   # 縦二軸のスケーリング

subset %>% 
  ggplot2::ggplot(ggplot2::aes(x = date)) + 
    ggplot2::geom_bar(ggplot2::aes(y = n), stat = "identity", width = 1.0,
                      fill = "dark gray", alpha = 1.0) + 
    ggplot2::geom_line(ggplot2::aes(y = ma7), linetype = "solid",
                       colour = "gray10", size = 0.5) + 
    # 第二軸を利用するグラフを描画する際はスケーリング調整する必要あり
    ggplot2::geom_line(ggplot2::aes(y = cum / sec_scale),
                       colour = "dark green", size = 1.0) +
    # 二軸表示のための軸属性の指定
    ggplot2::scale_y_continuous(
      # 第一軸のラベル(スケールは自動調整)
      name = "陽性者数(灰)・移動平均(黒)",
      # 第二軸の指定(第一軸にあわせたスケーリング)
      sec.axis = ggplot2::sec_axis(~ . * sec_scale,
                                   name = "累積陽性者数(濃緑)")) + 
    ggplot2::theme(
      axis.text.y.left = ggplot2::element_text(colour = "gray10"),
      axis.line.y.left = ggplot2::element_line(colour = "gray10"),
      axis.text.y.right = ggplot2::element_text(colour = "dark green"),
      axis.line.y.right = ggplot2::element_line(colour = "dark green")
    ) +
    # グラフ全体のラベル指定
    ggplot2::labs(title = title, subtitle = subtitle, caption = caption,
                  x = xlab, y = ylab)

 

同前日差

subset <- japan_daily %>% dplyr::mutate(y = diff)
title <- "【全国】前日差(単日)"
xlab <- ""
ylab <- "陽性者数"

subset %>% 
  ggplot2::ggplot(ggplot2::aes(x = date, y = y)) + 
    ggplot2::geom_line(colour = "dark green", alpha = 0.5) +
    ggplot2::labs(title = title, subtitle = subtitle, caption = caption,
                  x = xlab, y = ylab)

 

地方区分別集計

region_daily %>% 
  dplyr::filter(date == max(date))
subset <- region_daily %>% dplyr::mutate(key = region, y = n)
title <- "【地方別】陽性者数(単日)"
xlab <- ""
ylab <- "陽性者数"

subset %>% 
  ggplot2::ggplot(ggplot2::aes(x = date, y = y)) + 
    ggplot2::geom_bar(ggplot2::aes(fill = key), stat = "identity",
                      width = 1.0, alpha = 0.5) + 
    ggplot2::labs(title = title, subtitle = subtitle, caption = caption, 
                  x = xlab, y = ylab)

 

subset <- region_daily %>% dplyr::mutate(key = region, y = n)
title <- "【地方別】単日"
xlab <- ""
ylab <- "陽性者数"

subset %>% 
  ggplot2::ggplot(ggplot2::aes(x = date, y = y, colour = key)) + 
    ggplot2::geom_line(size = 0.5) +
    ggplot2::theme(legend.position = 'none') +
    ggrepel::geom_text_repel(ggplot2::aes(label = key),
                             data = subset(subset, date == max(date)),
                             nudge_x = 30, segment.alpha = 0.5, size = 4) + 
    ggplot2::lims(x = c(min(subset$date),
                        max(subset$date) + 45)) +
    ggplot2::labs(title = title, subtitle = subtitle, caption = caption,
                  x = xlab, y = ylab) 

 

同累計

subset <- region_daily %>% dplyr::mutate(key = region, y = cum)
title <- "【地方別】累計"
xlab <- ""
ylab <- "陽性者数"

subset %>% 
    ggplot2::ggplot(ggplot2::aes(x = date, y = y, colour = key)) + 
    ggplot2::geom_line(size = 1) +
    ggplot2::theme(legend.position = 'none') +
    ggrepel::geom_text_repel(ggplot2::aes(label = key),
                             data = subset(subset, date == max(date)),
                             nudge_x = 30, segment.alpha = 0.5, size = 4) + 
    ggplot2::lims(x = c(min(subset$date),
                        max(subset$date) + 45)) +
    ggplot2::labs(title = title, subtitle = subtitle, caption = caption,
                  x = xlab, y = ylab) 

 

同移動平均(7日)

subset <- region_daily %>% dplyr::mutate(key = region, y = ma7)
title <- "【地方別】移動平均(7日)"
xlab <- ""
ylab <- "陽性者数"

subset %>% 
    ggplot2::ggplot(ggplot2::aes(x = date, y = y, colour = key)) + 
    ggplot2::geom_line(size = 1) +
    ggplot2::theme(legend.position = 'none') +
    ggrepel::geom_text_repel(ggplot2::aes(label = key),
                             data = subset(subset, date == max(date)),
                             nudge_x = 30, segment.alpha = 0.5, size = 4) + 
    ggplot2::lims(x = c(min(subset$date),
                        max(subset$date) + 45)) +
    ggplot2::labs(title = title, subtitle = subtitle, caption = caption,
                  x = xlab, y = ylab) 

 

同移動平均(28日)

subset <- region_daily %>% dplyr::mutate(key = region, y = ma28)
title <- "【地方別】移動平均(28日)"
xlab <- ""
ylab <- "陽性者数"

subset %>% 
    ggplot2::ggplot(ggplot2::aes(x = date, y = y, colour = key)) + 
    ggplot2::geom_line(size = 1) +
    ggplot2::theme(legend.position = 'none') +
    ggrepel::geom_text_repel(ggplot2::aes(label = key),
                             data = subset(subset, date == max(date)),
                             nudge_x = 30, segment.alpha = 0.5, size = 4) + 
    ggplot2::lims(x = c(min(subset$date),
                        max(subset$date) + 45)) +
    ggplot2::labs(title = title, subtitle = subtitle, caption = caption,
                  x = xlab, y = ylab) 

 

同前日差

subset <- region_daily %>% dplyr::mutate(key = region, y = diff)
title <- "【地方別】前日差(単日)"
xlab <- ""
ylab <- "陽性者数"
ncol <- 3

subset %>% 
  ggplot2::ggplot(ggplot2::aes(x = date, y = y, colour = key)) + 
    ggplot2::geom_line(alpha = 0.75) +
    ggplot2::theme(legend.position = 'none') + 
    ggplot2::facet_wrap(~ key, ncol = ncol, scales = "free_y") + 
    ggplot2::labs(title = title, subtitle = subtitle, caption = caption,
                  x = xlab, y = ylab)

 

都道府県別

pref_daily %>% 
  dplyr::filter(date == max(date))
subset <- pref_daily %>% dplyr::mutate(key = pref, y = n)
title <- "【都道府県別】陽性者数(単日)"
xlab <- ""
ylab <- "陽性者数"

subset %>% 
    ggplot2::ggplot(ggplot2::aes(x = date, y = y, colour = key)) + 
    ggplot2::geom_line(size = 0.5) +
    ggplot2::theme(legend.position = 'none') +
    ggrepel::geom_text_repel(ggplot2::aes(label = key),
                             data = subset(subset, date == max(date)),
                             nudge_x = 30, segment.alpha = 0.5, size = 4) + 
    ggplot2::lims(x = c(min(subset$date),
                        max(subset$date) + 45)) +
    ggplot2::labs(title = title, subtitle = subtitle, caption = caption,
                  x = xlab, y = ylab) 

 

subset <- pref_daily %>% dplyr::mutate(key = pref)
title <- "【都道府県別】陽性者数(単日)"
xlab <- ""
ylab <- "陽性者数"
sec_scale <- 50
ncol <- 3

subset %>% 
  ggplot2::ggplot(ggplot2::aes(x = date)) + 
    ggplot2::geom_bar(ggplot2::aes(y = n, fill = key), stat = "identity",
                      alpha = 0.25, width = 1.0) + 
    ggplot2::geom_line(ggplot2::aes(y = ma7, colour = key),
                       linetype = "solid", size = 0.25) + 
    ggplot2::geom_line(ggplot2::aes(y = cum / sec_scale, colour = key)) +
    ggplot2::facet_wrap(~ key, ncol = ncol, scales = "free_y") + 
    ggplot2::theme(legend.position = 'none') + 
    ggplot2::scale_y_continuous(
      name = "陽性者数(棒)・移動平均(細線)",
      sec.axis = ggplot2::sec_axis(~ . * sec_scale,
                                    name = "陽性者数累計(太線)")) +
    ggplot2::labs(title = title, subtitle = subtitle, caption = caption,
                  x = xlab, y = ylab)

 

同前日差

subset <- pref_daily %>% dplyr::mutate(key = pref, y = diff)
title <- "【都道府県別】前日差(単日)"
xlab <- ""
ylab <- "陽性者数"
ncol <- 3

subset %>% 
  ggplot2::ggplot(ggplot2::aes(x = date, y = y, colour = key)) + 
    ggplot2::geom_line(alpha = 0.75) +
    ggplot2::theme(legend.position = 'none') + 
    ggplot2::facet_wrap(~ key, ncol = ncol, scales = "free_y") + 
    ggplot2::labs(title = title, subtitle = subtitle, caption = caption,
                  x = xlab, y = ylab)

 

年代別

subset <- ageBracket_daily %>% dplyr::mutate(key = ageBracket, y = n)
title <- "【年代別】陽性者数(単日)"
xlab <- ""
ylab <- "陽性者数"

subset %>% 
  ggplot2::ggplot(ggplot2::aes(x = date, y = y)) + 
    ggplot2::geom_bar(ggplot2::aes(fill = key), stat = "identity",
                      width = 1.0, alpha = 0.5) + 
    ggplot2::labs(title = title, subtitle = subtitle, caption = caption, 
                  x = xlab, y = ylab)

 

クラスター別

subset <- cluster_daily %>% dplyr::mutate(key = cluster, y = n)
title <- "【クラスター別】陽性者数(単日)"
xlab <- ""
ylab <- "陽性者数"


subset %>% 
  ggplot2::ggplot(ggplot2::aes(x = date, y = y)) + 
    ggplot2::geom_bar(ggplot2::aes(fill = key), stat = "identity",
                      width = 1.0, alpha = 0.5) + 
    ggplot2::labs(title = title, subtitle = subtitle, caption = caption, 
                  x = xlab, y = ylab)

 

対数化

地方区分(累計)

subset <- region_daily %>% dplyr::mutate(key = region, y = cum)
title <- "【地方別】累計"
xlab <- ""
ylab <- "陽性者数(常用対数)"

subset %>% 
    ggplot2::ggplot(ggplot2::aes(x = date, y = y, colour = key)) + 
    ggplot2::geom_line(size = 1) +
    ggplot2::theme(legend.position = 'none') +
    ggrepel::geom_text_repel(ggplot2::aes(label = key),
                             data = subset(subset, date == max(date)),
                             nudge_x = 30, segment.alpha = 0.5, size = 4) + 
    ggplot2::lims(x = c(min(subset$date),
                        max(subset$date) + 45)) +
    ggplot2::scale_y_log10() + 
    ggplot2::labs(title = title, subtitle = subtitle, caption = caption,
                  x = xlab, y = ylab) 

 

同移動平均(7日)

subset <- region_daily %>% dplyr::mutate(key = region, y = ma7)
title <- "【地方別】移動平均(7日)"
xlab <- ""
ylab <- "陽性者数(常用対数)"

subset %>% 
    ggplot2::ggplot(ggplot2::aes(x = date, y = y, colour = key)) + 
    ggplot2::geom_line(size = 1) +
    ggplot2::theme(legend.position = 'none') +
    ggrepel::geom_text_repel(ggplot2::aes(label = key),
                             data = subset(subset, date == max(date)),
                             nudge_x = 30, segment.alpha = 0.5, size = 4) + 
    ggplot2::lims(x = c(min(subset$date),
                        max(subset$date) + 45)) +
    ggplot2::scale_y_log10() + 
    ggplot2::labs(title = title, subtitle = subtitle, caption = caption,
                  x = xlab, y = ylab)

 

同移動平均(28日)

subset <- region_daily %>% dplyr::mutate(key = region, y = ma28)
title <- "【地方別】移動平均(28日)"
xlab <- ""
ylab <- "陽性者数(常用対数)"

subset %>% 
    ggplot2::ggplot(ggplot2::aes(x = date, y = y, colour = key)) + 
    ggplot2::geom_line(size = 1) +
    ggplot2::theme(legend.position = 'none') +
    ggrepel::geom_text_repel(ggplot2::aes(label = key),
                             data = subset(subset, date == max(date)),
                             nudge_x = 30, segment.alpha = 0.5, size = 4) + 
    ggplot2::lims(x = c(min(subset$date),
                        max(subset$date) + 45)) +
    ggplot2::scale_y_log10() + 
    ggplot2::labs(title = title, subtitle = subtitle, caption = caption,
                  x = xlab, y = ylab)

 


CC 4.0 BY-NC-SA, Sampo Suzuki